Contribution of Telomere Length to Systemic Sclerosis Onset: A Mendelian Randomization Study

Although previous studies have suggested a relationship between telomere shortening and systemic sclerosis (SSc), the association between these two traits remains poorly understood. The objective of this study was to assess the causal relationship between telomere length in leukocytes (LTL) and SSc using the two-sample Mendelian randomization approach, with the genome-wide association study data for both LTL and SSc. The results of inverse-variance weighted regression (OR = 0.716 [95% CI 0.528–0.970], p = 0.031) and the Mendelian randomization pleiotropy residual sum and outlier method (OR = 0.716 [95% CI 0.563–0.911], p = 0.035) indicate an association between telomere length and SSc. Specifically, longer genetically predicted LTL is associated with a reduced risk of SSc. Sensitivity tests highlight the significant roles of the variants rs10936599 and rs2736100 annotated to the TERC and TERT genes, respectively. Our findings suggest an influence of telomere length in leukocytes on the development of SSc.


Introduction
Systemic sclerosis (SSc) is an immune-mediated inflammatory disease (IMID) characterized by vascular damage, chronic inflammation and fibrotic involvement of connective tissues [1].SSc patients can be classified as limited cutaneous SSc (lcSSc) or diffuse cutaneous SSc (dcSSc) based on the extension of skin fibrosis [2].Similarly to other IMIDs, SSc is a complex disease involving interplay of genetic, epigenetic and environmental factors [2].
Telomeres are nucleoprotein structures located at chromosome ends that are implicated in the preservation of genome integrity and stability [10].These structures naturally undergo telomere shortening (TS), a process linked to inflammation and cellular senescence, both of which are implicated in the pathogenesis of SSc [2,[11][12][13][14].Of note, previous studies evaluating telomere length in leukocytes (LTL) in SSc have been undertaken; however, the findings have been inconsistent, and their significance regarding SSc pathogenesis remains unclear [15][16][17][18][19].
Mendelian randomization (MR) studies allow for exploring the potential causal relationship between a risk factor and a disease [20].These studies employ single-nucleotide polymorphisms (SNPs) as instrumental variables (IVs) and are based on the assumption that SNPs are associated with the risk factor but not with the disease or a confounding factor [20].In recent years, MR studies have been conducted to overcome the limitations of observational studies and to provide insights into the causal relationship between several IMIDs and multiple exposures [21].
In the present study, we aimed to investigate the causal relationship between LTL and SSc through an MR study.

Results
Seven IVs were used in our study based on their association with LTL [22].This information and their respective effect estimators for the SSc and SSc clinical subtypes are shown in Table 1 and Table S1, respectively.The selected IVs were strong enough to avoid weak bias, as indicated by the F-statistic value of 67.03 [23].Since none of the heterogeneity tests were significant, the inverse-variance weighted (IVW) fixed-effects (FE) method was used for the three datasets analyzed (Tables 2 and 3).In addition, MR-Egger showed no evidence of horizontal pleiotropy in any of the datasets, and the MR pleiotropy residual sum and outlier (MR-PRESSO) method detected no outliers (Tables 2 and 3).Our results for IVW-FE show an association between genetically predicted longer LTL and a reduced risk of SSc (OR = 0.716 [95% CI 0.528-0.970],p = 0.031).In addition, three other MR methods, maximum likelihood, weighted median and MR-PRESSO, showed the consistency of our results, with a significant association and the same direction of effect (Table 2).Interestingly, even when weighted mode regression did not reach statistical significance, the estimate was in the same direction as the other methods.MR-Egger regres-sion results were not statistically significant (OR = 0.397 [95% CI 0.094-1.670],p = 0.263).Table 2 summarizes the MR results for LTL and SSc, indicating a negative relation between genetically predicted telomere length and SSc.Regarding leave-one-out (LOO) sensitivity analysis, we observed that individually removing two SNPs (rs10936599 and rs2736100) from the analysis resulted in the loss of the significance, maintaining the direction of the effect (p = 0.241 and p = 0.102, respectively; Table S2).Lastly, we investigated the specific associations for LTL with the two main clinical subtypes of the disease.Our results show a significant association between genetically predicted longer LTL and a decreased risk of lcSSc for IVW-FE with an OR of 0.669 [95% CI 0.468-0.956]and a p-value of 0.027.In addition, the maximum likelihood, weighted median and MR-PRESSO methods were also statistically significant (Table 3).In contrast, the results in dcSSc showed no statistically significant association with LTL for any of the MR methods applied (IVW-FE: OR = 0.771 [95% CI 0.481-1.237],p = 0.281; Table 3).Noteworthy, LOO sensitivity analysis results for lcSSc revealed the loss of significance for the same two SNPs as observed for the complete SSc dataset, rs10936599 p = 0.234 and rs2736100 p = 0.081 (Table S3).

Discussion
In the current study, we investigated the possible causal relationship of LTL with the risk of SSc, and its main clinical subtypes, using the MR methodology and the largest SSc GWAS in Europeans [3].Our primary findings show a significant association between genetically predicted longer LTL and a reduced risk of SSc, shedding light on the direction of influence between these two traits.
Previous observational studies in LTL in SSc reached conflicting results [15][16][17].However, the largest and most recent study described a higher proportion of patients with SSc with shorter telomeres compared to controls [15].Consistent with this observation, our data indicate that genetically predicted longer LTLs are associated with a reduced risk of SSc (IVW-FE: OR = 0.716 [95% CI 0.528-0.970],p = 0.031; Table 2).Furthermore, owing to the methodology applied, our study provides evidence that LTL influences the risk of developing SSc rather than the disease progression affecting LTL.These results align with the data reported in a MR study for LTL performed for an IMID primarily affecting the skin, psoriasis [24], and another connective tissue disease, rheumatoid arthritis [25].However, our results significantly contrast with the results of a similar study in systemic lupus erythematosus [26].These discrepancies may be explained by the high complexity of connective tissue disorders.
Regarding the relationship of LTL with the main clinical subtypes, we observed a significant association between LTL and lcSSc, whereas the analysis for dcSSc did not reach statistical significance (Table 3).This observation may be attributed to variations in statistical power between clinical subtypes, as shown in Table S4.Several studies have analyzed LTL in lcSSc and dcSSc, yielding conflicting results [15,18,19] that may be related to the different sample sizes, normal telomere length (TL) reference or methods used to measure TL.These inconsistencies underscore the necessity for further and larger studies to understand LTL in lcSSc and dcSSc.
In examining the sensitivity analyses, we found that two SNPs (rs10936599 and rs2736100) were the primary contributors to the obtained results.Interestingly, these SNPs and their corresponding annotated genes, TERC and TERT, were associated with idiopathic forms of interstitial lung disease [27,28].Moreover, it is worth nothing that telomere shortening is an established risk factor for idiopathic pulmonary fibrosis [14,24] and has also been associated with different pulmonary features in SSc patients [15,18].Taking the above into consideration, the contribution of LTL to SSc seems to be especially relevant in the severe pulmonary affection of the disease.Unfortunately, the clinical data regarding pulmonary affection in our large GWAS cohort are limited, and we could not evaluate this in the present study.
Our results underscore the potential involvement of LTL in SSc development.However, the mechanism by which LTL affects SSc remains unknown.In this context, we hypothesize some potential pathways.First, short TL could contribute to the dysregulation of the immune system by affecting T-cell numbers and receptor diversity [28].Notably, CD28-negative T cells that have a proinflammatory profile with alternative receptors, cytolytic properties and shorter TL [13] have been reported to be increased in the blood and skin of SSc patients [29].Additionally, several studies have related short TL and telomere dysfunction with increased proinflammatory cytokines [30][31][32], which could contribute to the immune imbalance of the disease.Furthermore, the presence of autoantibodies against a telomere-related protein (TERF1) in some SSc patients and its association with shorter LTL [33] highlight the interplay between telomeres and immune response.

Materials and Methods
To determine the possible causal relationship between LTL and SSc and its main clinical subtypes, lcSSc and dcSSc, we carried out a two-sample MR (2SMR) study (Figure 1).This approach allows the use of two non-overlapping datasets, one for the exposure and another for the outcome, to determine the causal relationship between them [34].

Genetic Data Sources and IV Selection
Genetic association data for seven non-palindromic independent SNPs, associated at genome-wide significance level (p < 5 × 10 −8 ) with LTL, were obtained from a GWAS meta-analysis comprising 37,684 European-descent individuals [22].Summary statistics from the largest SSc GWAS in Europeans (9095 patients with SSc and 17,584 controls) [3] were used to extract the association estimates for the outcome.A summary of these SNPs as IVs and their size effect on both LTL and SSc can be found in Table 1.In addition, association data specific for lcSSc and dcSSc were also used (Table S1).The sample size of the cohort and the clinical subtypes are shown in Table S5.
To ensure the adequacy of the IVs, Phenoscanner [35,36] and LDtrait [37,38] were used to verify the absence of a direct association of the IVs with the outcome.The strength of the selected IVs was evaluated using the F statistic, with an F value greater than 10 indicating sufficient strength to avoid weak bias [23].Furthermore, we calculated the statistical power of our analysis following the methodology of Brion et al. [39].

Genetic Data Sources and IV Selection
Genetic association data for seven non-palindromic independent SNPs, associated at genome-wide significance level (p < 5 × 10 −8 ) with LTL, were obtained from a GWAS metaanalysis comprising 37,684 European-descent individuals [22].Summary statistics from the largest SSc GWAS in Europeans (9095 patients with SSc and 17,584 controls) [3] were used to extract the association estimates for the outcome.A summary of these SNPs as IVs and their size effect on both LTL and SSc can be found in Table 1.In addition, association data specific for lcSSc and dcSSc were also used (Table S1).The sample size of the cohort and the clinical subtypes are shown in Table S5.
To ensure the adequacy of the IVs, Phenoscanner [35,36] and LDtrait [37,38] were used to verify the absence of a direct association of the IVs with the outcome.The strength of the selected IVs was evaluated using the F statistic, with an F value greater than 10 indicating sufficient strength to avoid weak bias [23].Furthermore, we calculated the statistical power of our analysis following the methodology of Brion et al. [39].

MR Analysis
MR analyses were performed using the 2SMR approach with the R package "TwoSampleMR" [34], and a significant association was determined at a p < 0.05.The IVW, the maximum likelihood, the MR-Egger and the MR-PRESSO methods were selected for the analysis.
The IVW method combines the effects of all the IVs and, by assuming the validity or invalidity of all the SNPs, sets the global pleiotropy to zero [34].The selection between IVW FE or random effects was determined by the p-value of the heterogeneity test, where a p > 0.05 indicated the use of the fixed-effects model.The maximum likelihood approach is based on the direct maximization of the likelihood and provides more accurate confidence intervals than the IVW when there is some uncertainty in the genetic associations [40].We also employed MR-Egger regression, as it has the capacity to estimate causality using weak or invalid IVs and to provide an estimation of the pleiotropy [41].MR-PRESSO

MR Analysis
MR analyses were performed using the 2SMR approach with the R package "TwoSam-pleMR" [34], and a significant association was determined at a p < 0.05.The IVW, the maximum likelihood, the MR-Egger and the MR-PRESSO methods were selected for the analysis.
The IVW method combines the effects of all the IVs and, by assuming the validity or invalidity of all the SNPs, sets the global pleiotropy to zero [34].The selection between IVW FE or random effects was determined by the p-value of the heterogeneity test, where a p > 0.05 indicated the use of the fixed-effects model.The maximum likelihood approach is based on the direct maximization of the likelihood and provides more accurate confidence intervals than the IVW when there is some uncertainty in the genetic associations [40].We also employed MR-Egger regression, as it has the capacity to estimate causality using weak or invalid IVs and to provide an estimation of the pleiotropy [41].MR-PRESSO was implemented due to its capability to detect outliers and to provide an outlier-free estimate [42].Additionally, weighted median and weighted mode were applied to obtain an unbiased estimate in the presence of some invalid IVs, weighting the contribution of the IVs [43,44].The median method requires at least half of IVs to be valid, while the mode method assumes the validity of IVs within the larger group of IVs based on their similarity [43,44].
In order to evaluate the effect of each SNP, we carried out a LOO sensitivity analysis, which provides the IVW estimate sequentially excluding one of the IVs at a time [34], applying the tool implemented in the "TwoSampleMR" package (https://elifesciences.org/ articles/34408), accessed on 26 September 2023 [34].

Conclusions
In conclusion, our study points towards a possible causal relationship between increased LTL and a reduced risk of SSc.It is important to interpret these results with caution until further investigation is conducted in this area, acknowledging the inherent limitations of MR studies in fully assessing the assumptions related to the exclusive impact of the IVs on the outcome through the risk factor.Future studies with the aim of clarifying specific associations and evaluating possible mechanisms are warranted.

Supplementary Materials:
The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ijms242115589/s1.Institutional Review Board Statement: An ethical protocol was prepared, consensus was reached across all partners, academic and industrial, translated into all participants' languages, and approved by each of the local ethical committees of the clinical recruitment centers.The studies adhered to the standards set by the International Conference on Harmonization and Good Clinical Practice (ICH-GCP) and to the ethical principles that have their origin in the Declaration of Helsinki (2013).The protection of the confidentiality of records that could identify the included subjects is ensured as defined by the EU Directive 2001/20/EC and the applicable national and international requirements relating to data protection in each participating country.For the systemic sclerosis data, as previously described (López-Isac E et al.GWAS for systemic sclerosis identifies multiple risk loci and highlights fibrotic and vasculopathy pathways.Nat Commun.2019.https://doi.org/10.1038/s41467-019-12760-y).CSIC's ethics committee approved the study protocol, and written informed consent was obtained in accordance with the tenets of the Declaration of Helsinki.In the case of the leukocyte telomere length, GWAS summary statistics data were public at the time of the study.
Informed Consent Statement: Not applicable.

Figure 1 .
Figure 1.Framework of the 2SMR study for leukocyte telomere length and systemic sclerosis.The 2SMR approach utilizes data from two independent GWAS and employs SNPs as instrumental variables.The selected SNPs are associated with the exposure, influence the outcome only through the exposure and should not be associated with any confounder.Parts of the figure were drawn by using and modifying pictures from Servier Medical Art.Servier Medical Art by Servier is licensed under a Creative Commons Attribution 3.0 Unported License (https://creativecommons.org/licenses/by/3.0/),accessed on 26 September 2023.2SMR: two-sample Mendelian randomization; GWAS: genome-wide association studies; IV: instrumental variable; SNP: single-nucleotide polymorphism.

Figure 1 .
Figure 1.Framework of the 2SMR study for leukocyte telomere length and systemic sclerosis.The 2SMR approach utilizes data from two independent GWAS and employs SNPs as instrumental variables.The selected SNPs are associated with the exposure, influence the outcome only through the exposure and should not be associated with any confounder.Parts of the figure were drawn by using and modifying pictures from Servier Medical Art.Servier Medical Art by Servier is licensed under a Creative Commons Attribution 3.0 Unported License (https://creativecommons.org/licenses/ by/3.0/),accessed on 26 September 2023.2SMR: two-sample Mendelian randomization; GWAS: genome-wide association studies; IV: instrumental variable; SNP: single-nucleotide polymorphism.

Author Contributions:
Study design: J.M. and M.A.-H.Clinical interpretation of the data: A.G.-D.-C., N.O.-C., J.L.C. and C.P.S.-A.Data analysis and interpretation: I.R.-M.and G.V.-M.Writing-original draft preparation: I.R.-M., G.V.-M., J.M. and M.A.-H.All authors have read and agreed to the published version of the manuscript.Funding: This work was supported by Redes de Investigación Cooperativa Orientadas a Resultados en Salud (RICORS) (RD21/0002/0039) funded by Instituto de Salud Carlos III (ISCIII) from the Spanish Ministry of Science and Innovation.M.A.H. is a recipient of a Miguel Servet fellowship (CP21/00132) from ISCIII.G.V.M was funded by Grant PRE2019-087586 funded by MCIN/AEI/10.13039/501100011033and by "ESF Investing in your future".I.R.M was supported by the program FPU (FPU21/02746) funded by the Spanish Ministry of University.

Table 1 .
Genetic associations of the selected instrumental variables with LTL and SSc.

Table 3 .
Association between genetically predicted LTL and risk of SSc main clinical subtypes.